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Abstract 

A review of the coupled cluster method (CCM) applied to lattice quantum spin systems is 
presented here. The CCM formalism is explained and an application to the spin-half XXZ model 
on the square lattice is presented. Low orders of approximation are carried out analytically and 
a high-order CCM formulation is presented. Results for the SUB2 approximation are carried out 
numerically for SUB2-m and analytically for SUB2. It is found that SUB2-m results converge 
rapidly to the full SUB2 solution, including new results for the SUB2-m limiting points compared 
to the SUB2 critical point. Results for the ground-state energy and the sublattice magnetisation 
are presented. A study of the excitation spectrum of this model at the SUB2 critical point is 
given. Indeed, the shape of the excitation spectrum at the SUB2 critical point is identical to that 
predicted by spin-wave theory for the isotropic model, albeit with a multiplicative factor of 1.1672. 
This result compares very well to results of cumulant series expansions and Monte Carlo simulation 
that again predict a similar shape for the excitation spectrum, but with multiplicative factors of 
1.18 and 1.21ib0.03, respectively. Results for the isotropic Heisenberg model on the square lattice 
for the spin-one antiferromagnet and the spin-one/spin-half ferrimagnet are also given. 
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I. INTRODUCTION 



Lattice quantum spin systems at zero temperature may be treated using a variety of 



methods, including exact solutions 



mm 



in certain the number of lattice sites, 



A^, goes to N ^ oc. However, the majority of cases for quantum spin systems may not be 
solved by using such exact techniques and so we must use approximate methods. An advance 
in the treatment of low- dimensional quantum systems has been seen with the introduction 
of the density matrix renormalisation group (DMRG) technique jSj. This method provides 
accurate results consistently for a variety of one-dimensional or quasi-one-dimensional quan- 
tum systems. However, it is still perhaps fair to say that the power of the DMRG method 
still remains largely untapped for systems of two or higher spatial dimensions, although we 



note that recent advances are encouraging |6| . The properties of "unfrustrated" lattice quan- 
tum spin systems of two spati_al_dimensions may however be considered using the quantum 
Monte Carlo (QMC) method 



m 



We note that QMC calculations offer the possibility of a 
level of accuracy limited only by the computing power available. We note that errors scale in 
a statistically well-understood manner with the length of the MC simulation. However, the 
QMC method is severely limited by the infamous sign problem, which at zero temperature is 
a reflection of quantum "frustration" in the spin system. There are many other approximate 
techniques which may be applied, such as spin- wave theory j^, exact diagonalisations jl^ . 
and cumulant series expansions although these techniques are not discussed here. 

We wish to focus here on another approximate technique called the coupled cluster 
method (CCM). The CCM often provides accurate results even in the presence of very 
strong frustration. Note that the CCM technique 
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20| is a 



well-known method of quantum many-body theory (QMBT). We note that the CCM has 



been applied with much success over t 



at zero ternperature 



mm, 

mm 
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In particular, the use of computer-algebraic implementations 
44j | of the CCM for quantum systems of large or infinite numbers of particles has 
largely been found to be very effective with respect to these spin-lattice problems. 

We present a brief description of the CCM formalism. We then describe the application of 
the method to the spin-half XXZ model for the square lattice at zero temperature using two 
model states. The practicalities of using the CCM are illustrated by considering the LSUB2 
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and SUB2 approximations in detail, and their analytical solution is shown. The method of 
carrying out high-order CCM calculations using computational approaches for the localised 
LSUBm and SUBm-m approximation schemes is described in detail. A calculation for the 
excitation spectra of the spin-half XXZ model using the SUB2 approximation in the ground 
state is presented. Results for excitation energies determined using localised approximation 
schemes to high order is also presented. We then broaden our treatment of the Heisenberg 
model on the square lattice to include quantum spin numbers of s > 1/2. The results are 
found to be in good agreement with those results of the best of other approximate methods, 
where they exist. We conclude with a discussion of the implications of these results and of 
future applications of the CCM. 



II. THE CCM FORMALISM 



A brief description of the normal coupled cluster method (NCCM) formalism is now 



provided, alth oug h the interested reader is referrec 
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for further details. The exact ket and bra ground-state energy eigenvectors, 



and of a general many-body system described by a Hamiltonian H 



are parametrised within the single-reference CCM as follows: 



(1) 



1^) = 6^1$) 

(§1 = ($|5e-^ 



s = Y.SiCt, 



(2) 



The single model or reference state |$) is required to have the property of being a cyclic 
vector with respect to two well-defined Abelian subalgebras of multi-configurational creation 
operators {Cf} and their Hermitian-adjoint destruction counterparts {Cj = {C^Y}. Thus, 
1$) plays the role of a vacuum state with respect to a suitable set of (mutually commuting) 
many-body creation operators {C/}. Note that C7|$) = 0, V / 7^ 0, and that Cq = 1, 
the identity operator. These operators are furthermore complete in the many-body Hilbert 
(or Fock) space. Also, the correlation operator S is decomposed entirely in terms of these 
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creation operators {Cf}, which, when acting on the model state ({C/|$)}), create exci- 
tations from it. We note that ahhough the manifest Hermiticity, = is 
lost, the normalisation conditions (^^I^E') = (*!&|^) = ('^'1*^') = 1 are explicitly imposed. The 
correlation coefficients {Sj, Sj} are regarded as being independent variables, and the full set 
{Sj, Sj} thus provides a complete description of the ground state. For instance, an arbitrary 
operator A will have a ground-state expectation value given as, 

A = = (<l>|5e-^Ae^|$) = A {{Si, Si}) . (3) 

We note that the exponentiated form of the ground-state CCM parametrisation of Eq. Q 
ensures the correct counting of the independent and excited correlated many-body clusters 
with respect to |$) which are present in the exact ground state It also ensures the 

exact incorporation of the Goldstone linked-cluster theorem, which itself guarantees the 
size-extensivity of all relevant extensive physical quantities. We also note that any operator 
in a similarity transform may be written as 

A ^ e-'Ae' = A+[A,S] + ^[[A, S],S] + --- (4) 

The determination of the correlation coefficients {Si, Si} is achieved by taking appropriate 
projections onto the ground-state Schrodinger equations of Eq. (^. Equivalently, they may 
be determined variationally by requiring the ground-state energy expectation functional 
H{{Si, Si}), defined as in Eq. Q, to be stationary with respect to variations in each of the 
(independent) variables of the full set. We thereby easily derive the following coupled set of 
equations, 

6H/6S1 = ^ {<l>\Cfe-^He^\^) = 0, V / ^ ; (5) 
6H/6S1 = ^ {<!>\Se"^[H, C/]e^|<l>) = 0, V / ^ . (6) 

Equation (jSJ also shows that the ground-state energy at the stationary point has the simple 
form 

E, = E,{{Si}) = ($|e- . (7) 

It is important to realize that this (bi-)variational formulation does not lead to an upper 
bound for Eg when the summations for S and S in Eq. Q are truncated, due to the 
lack of exact Hermiticity when such approximations are made. However, one can prove 
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that the important Hellmann-Feynman theorem is preserved in all such approximations. 
We note that Eq. (0) represents a coupled set of non-linear multinomial equations for the 
c-number correlation coefficients {Sj}. The nested commutator expansion of the similarity- 
transformed Hamiltonian 

H ^ e-'He' = H+[H,S] + ^[[H, SiS] + --- (8) 

and the fact that all of the individual components of S in the sum in Eq. ^ commute 
with one another, together imply that each element of S in Eq. Q is linked directly to the 
Hamiltonian in each of the terms in Eq. Thus, each of the coupled equations (0) is of 
linked cluster type. Furthermore, each of these equations is of finite length when expanded, 
since the otherwise infinite series of Eq. (jH} will always terminate at a finite order, provided 
only (as is usually the case) that each term in the second-quantised form of the Hamiltonian, 
if, contains a finite number of single-body destruction operators, defined with respect to 
the reference (vacuum) state |$). Hence, the CCM parametrisation naturally leads to a 
workable scheme which can be efficiently implemented computationally. It is important to 
note that at the heart of the CCM lies a similarity transformation, in contrast with the 
unitary transformation in a standard variational formulation in which the bra state is 
simply taken as the explicit Hermitian conjugate of In the case of spin-lattice problems 
of the type considered here, the operators Cf become products of spin-raising operators 
s~l over a set of sites {fc}, with respect to a model state \(^) in which all spins points 
"downward" in some suitably chosen local spin axes. The CCM formalism is exact in the 
limit of inclusion of all possible such multi-spin cluster correlations for 5* and 5*, although 
in any real application this is usually impossible to achieve. It is therefore necessary to 
utilise various approximation schemes within 5* and S. The three most commonly employed 
schemes previously utilised have been: (1) the SUBn scheme, in which all correlations 
involving only n or fewer spins are retained, but no further restriction is made concerning 
their spatial separation on the lattice; (2) the SUBra-m sub-approximation, in which all 
SUBn correlations spanning a range of no more than m adjacent lattice sites are retained; 
and (3) the localised LSUBm scheme, in which all multi-spin correlations over all distinct 
locales on the lattice defined by m or fewer contiguous sites are retained. The problem 
of solving for these types of approximation schemes using analytical and computational 
approaches is discussed below. 
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An excited-state wave function, |\I/e), is determined by linearly applying an excitation 
operator to the ket-state wave function of Eq. (j21), such that 

|^e)=X^e^|$) . (9) 

This equation may now be used to determine the low-lying excitation energies, such that the 
Schrodinger equation, Hl'^e) = -f'el^e), may be combined with its ground-state counterpart 
of Eq. ([H) to give the result, 

eeX"|<l>) = e-'^[/7,X"]e^|$) , (10) 

where e^. = Ef. — Eg is the excitation energy. By analogy with the ground-state formalism, 
the excited-state correlation operator is written as, 

x^ = Y.x!Ct , (11) 

where the set {C/} of multi-spin creation operators may differ from those used in the 
ground-state parametrisation in Eq. Q if the excited state has different quantum numbers 
than the ground state. We note that Eq. (jllj) implies the overlap relation ($|\l/e) = 0. By 
applying ($|C7 to Eq. (fTUI) we find that, 

e^Xf = {^Cie-^[H,X^]e^\^) ,V/^0 , (12) 

which is a generalised set of eigenvalue equations with eigenvalues eg and corresponding 
eigenvectors Xf, for each of the excited states which satisfy (•I'l^e) = 0. We note that 
lower orders of approximation may be determined analytically for both the ground and 
excited states. Examples of applying the LSUB2 and SUB2 approximations to the spin-half 
square-lattice XXZ model are given later in order to show clearly how this is performed. 
However, it rapidly becomes clear that analytical determination of the CCM equations for 
higher orders of approximation is impractical. We therefore employ computer algebraic 
techniques in order efficiently to determine and solve the CCM ket- and bra-state equations 
(discussed below). The bra-state equations may be determined easily thereafter and the 
ket- and bra-state equations are readily solved using standard techniques for the solution 
of coupled polynomial equations (e.g., the Newton- Raphson method). The excited-state 
eigenvalue equations may be also determined in an analogous manner, and, although this 
is not strictly necessary, we restrict the level of approximation to the same for the excited 
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state as for the ground state in calculations presented here. A full exposition of the details 



in a pply in g the CCM to high orders of approximation is given for the ground state in Refs 



31 



36 



44| and for excited states in Ref. j39|. Note that the results of SUBm-m and 
LSUBm approximation schemes may be extrapolated to the exact limit, m oo, using 
various "heuristic' approaches. How to do this is not discussed here, although the interested 



reader is referred to Refs. 



36 



44| for more details. 



III. THE XXZ MODEL ON THE SQUARE LATTICE 

We wish to apply the CCM to the spin-half XXZ model on the square lattice in order 
to illustrate how one applies the CCM to a practical problem. We note that this system is 
unfrustrated and the XXZ Hamiltonian is specified as follows, 

^ = E[^?^j+^r4+M^i] ' (13) 

where the sum on counts all nearest-neighbour pairs once. The Neel state is the ground 
state in the trivial Ising limit A ^ oo, and a phase transition occurs at (or near to) A = 1. 
Indeed, the ground state demonstrates Neel-like order in the 2;-direction for A > 1 and a 
similar x-y planar phase for —1 < A < 1. The system is ferromagnetic for A < —1. We 
note that approximately 61% of the classical ordering remains in the quantum Heisenberg 
model at A = 1. 



A. The CCM applied to the XXZ model 

We turn now to the choice of model state |$) and the operators {Cf} for the case of spin- 
half quantum antiferromagnets on bipartite lattices. In the regime where the corresponding 
classical limit is described by a Neel-like order in which all spins on each sublattice are 
separately aligned in some global spin axes, it is convenient to introduce a different local 
quantisation axis and different spin coordinates on each sublattice. This is achieved by a 
suitable rotation in spin space, so that the corresponding reference state becomes a fully 
aligned ("ferromagnetic") state, with all spins pointing along, say, the negative z-axis in 
the corresponding local axes. Such rotations are cannonical tranformations that leave the 
spin commutation relations unchanged. In the same local axes, the configuration indices 
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J — > {fci, ^2, • • •, ^m}, a set of site indices, such that Cf ■^ki-^k2 ' ' ' -^kM^ where = sl±isl 
are the usual spin-raising and spin-lowering operators at site k. For the Hamiltonian of Eq. 
(|13p we first choose the 2;-aligned Neel state as our reference state (which is the exact ground 
state for A ^ oo, and is expected to be a good starting point for all A > 1, down to the 
expected phase transition at A = 1). We perform a rotation of the up-pointing spins by 
180° about the y-axis, such that —s^, s^, — > —s^ on this sublattice. The 

Hamiltonian of Eq. (fT^ may thus be written in these local coordinates as, 

H = -\Y.i4s^ + s■s~ + 2^sls^^, . (14) 

There is never a unique choice of model state |$). Indeed, our choice should be guided by any 
physical insight available to us concerning the system or, more specifically, that particular 
phase of it which is under consideration. In the absence of any other insight into the quantum 
many-body system, we may sometimes be guided by the behaviour of the corresponding 
classical system. The XXZ model under consideration provides just such an illustrative 
example. Thus, for A > 1 the classical Hamiltonian of Eq. ()13|1 on the 2D square lattice 
(and, indeed, on any bipartite lattice) is minimized by a perfectly antiferromagnetically 
Neel-ordered state in the 2;-direction, and we have already utilised this information in the 
preceding subsections. However, the classical ground-state energy is minimized by a Neel- 
ordered state with spins pointing along any direction in the xy plane, say along the x-axis 
for — 1 < A < 1. Thus, in order to provide CCM results in the region — 1 < A < 1, we now 
take this state to be our model state and we shall refer to it as the "planar" model state. 

In order to produce another "ferromagnetic" model state for the planar model state in the 
local frames, we rotate the axes of the left-pointing spins (i.e., those pointing in the negative 
x-direction) in the planar state by 90° about the y-axis, and the axes of the corresponding 
right-pointing spins by —90° about the y-axis. (Note that the positive 2;-axis is defined 
here to point upwards and the positive a;-axis is defined to point rightwards.) Thus, the 
transformations of the local axes are described by 

s"^s" , s^^s^ , s^^-s" (15) 

for the left-pointing spins, and by 

, s^^s^ , s^^s" (16) 
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for the right-pointing spins. The transformed Hamiltonian of Eq. ()13p may now be written 
in these local axes as 

^ = -7 E [(A + ^)i44 + S-S-) + (A - l){sts- + s-sf) + 4sts^] , (17) 

We note that all of the CCM correlation coefficients are zero at A = — 1 because the model 
state is an exact ground eigenstate of the Hamiltonian of Eq. (jl7|) at this point. Hence, we 
track the CCM solution for the planar model state from this "trivial" point at A = —1. 

The results presented below are based on the SUB2 approximation scheme and the 
localised SUBm-m and LSUBm schemes. We include all fundamental configurations, 
I {ki,k2,- ■ -kn}, with n < m, which are distinct under the point and space group 
symmetries of both the lattice and the Hamiltonian. The numbers, Np and Np^, of such 
fundamental configurations for the ground and excited states, respectively, may be further 
restricted by the use of additional conservation laws. For example, the Hamiltonian of Eq. 
(fT^ commutes with the total uniform magnetisation, = J2k where the sum on k runs 
over all lattice sites. The ground state is known to lie in the = subspace, and hence we 
exclude configurations with an odd number of spins or with unequal numbers of spins on the 
two equivalent sublattices for the ^-aligned model state. A similar condition is imposed on 
clusters for the ground state based on the planar model state. Similarly for the excited states 
with the z-aligned model state, since we are only interested in the lowest-lying excitation, 
we restrict the choice of configurations to those with = ±1. 

B. The LSUB2 approximation for the spin-half, square-lattice XXZ model for the 
z-aligned model state 

We start the LSUB2 calculation by specifying the commutation relations [sf , s|] = 
Ts't^i,k and [s^, s^] = 2sl5i^k- We again note that the similarity transform may be expanded 
as a series of nested commutators in Eq. We write the LSUB2 ket-state operator in the 
following simple form for the spin-half linear chain model, 

h ^ 

S=iEE4st,, (18) 

where i runs over all sites on the square lattice and p runs over all nearest-neighbour lattice 
vectors. Note that hi is the sole ket-state correlation coefficient. In this approximation 



we may therefore determine similarity transformed versions of the spin operators exphcitly, 
given by 



St = st + hY.4st,, (19) 



Si = si - 2bi J2 sf^ip -blJ2 44+pA+P2 

P Pl,P2 

We note that the otherwise infinite series of operators in the expansion of the similarity 
transform terminates to finite order. We also note that = for any lattice site 

(which is true only for spin-half systems), and this is implicitly assumed in the last of Eqs. 
()19|) . Clearly we may also write the similarity transformed version of the Hamiltonian as 

^ = E [St St + S-SJ + 2AStr,] . (20) 

We may now substitute the expressions for the spin operators in Eq. ()19|) into the above 
expression. The ground-state energy is given by 

f = -^{A + 26i} . (21) 

We note that our expression for the ground-state energy is size-extensive (i.e., it scales 
linearly with A^), as required by the Goldstone theorem which is obeyed by the NCCM. 
Furthermore, this expression terminates to finite order, as for the similarity transformed 
versions of spin operators. Finally, we note that any other non-trivial choice for S will 
always yield this expression for the ground-state energy. The task is now to find bi and we 
note that if we could include all possible spin correlations in S then we would obtain an exact 
result for the ground-state energy. However, this is found to be impossible to achieve for 
most cases in practice, and we make an approximation (such as the LSUB2 approximation 
presented here). The LSUB2 ket-state equation is given by 

5bl + 6A6i -1=0, (22) 

which therefore implies that the LSUB2 ground-state energy maybe written explicitly in 
terms of A as 
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We note that this expression gives the correct result in the Ising hmit A — > oo. We again 
note that the bra state does not manifestly have to be the Hermitian conjugate of the ket 
state, and we note that the bra-state correlation operator for the LSUB2 approximation is 
given by, 

^=1 + |EEw. ' (24) 

where the index j runs over all sites on the linear chain and bi is the sole bra-state correlation 
coefficient in the LSUB2 approximation. In order to determine the bra-state equation, we 
now explicitly determine H{{Si,Si}), 

H = - — {A + 2bi) + Nh(6Abi + 5bl-l) , (25) 

such that LSUB2 bra-state equation is given from dH/dbi = as 

6A6i + IO6161 -1 = 0, (26) 

which gives bi — |(9A^ -|- 5)"^/^. Finally, we note that once the values for the bra- and ket- 
state correlation coefficients have been determined (at a given level of approximation) then 
we may also obtain the values for expectation values, such as the sublattice magnetisation 
given by 

M^-^J^\j:sm--^{^\Se-'{f:st)e'm . (27) 

The sublattice magnetisation is written here in terms of the "rotated" spin coordinates. We 
note that this is given by, 

Mlsub2 = 1 - 86161 , 

= 1[1+ ] (28) 

for the LSUB2 approximation. 

C. The SUB2 approximation for the spin- half, square- lattice XXZ model of the 
z-aligned model state 

The SUB2 approximation allows us to include all possible two-spin correlations in our 
wave function. We note that the SUB2 ket-state operator is given by 

1 ^ 

'^=9EEM+4+. , (29) 

11 



and that the index i runs over all sites on the linear chain. Furthermore, the index r runs over 
all lattice vectors which connect one sublattice to the other and bj. is its corresponding SUB2 
ket-state correlation coefficient for this vector. We again determine a similarity transformed 
version of the spin operators and we are able to determine the SUB2 equations, given by 

J2{{l + 2Abi + 2bl)6p,r-2{A + 2bi)br + J2br'+p+pA-r'-pi} =0 , (30) 
P r' 

where p runs over all (2D) nearest-neighbour lattice vectors and pi is any one of these lattice 
vectors. Equation pUj) may now be solved by employing a sublattice Fourier transform, given 
by 

r(g)=^e'^-^6. , (31) 

r 

where r again is a lattice vector (i.e., + Vy is an odd integer number for the 2D square 
lattice) which connects the different sublattices. This expression has an inverse given by 

br = ^ dq^ dqyCos{r^q^)cos{ryqy)T{q) . (32) 
vr^ JO Jo 

The SUB2 equations (jHUjl and Eq. (j!TT|) therefore lead to an expression for T{q) given by 

r(g) = :^[i±Vi-fcV(9)] , (33) 

where K = A + 2bi, = {1 + 2A6i + 2b1)/K'^, and 7(g) = 0.5(cosqx + cosqy). (Note 
that we choose the negative solution in Eq. ()33|) in order to reproduce results in the trivial 
limit A ^ 00.) These equations now yield a self-consistency requirement on the variable 
bi and they may be solved iteratively at a given value of A. Indeed, we know that all 
correlation coefficients must tend to zero (namely, for SUB2: br ^ 0, V r) as A ^ cxd and 
we track this solution for large A by reducing A in small successive steps. We find that the 
discriminant in Eq. ()33j) becomes negative at a critical point, Ac ~ 0.7985. This is a strong 
indication that the CCM critical point is detecting the known quantum phase transition in 
the system at A = 1. Furthermore, the SUB2 approximation for the ground state may be 
used in conjunction with a SUBl approximation for the excited state operator in Eq. 
(jlip in order to determine the excitation spectrum. We note that the excitation spectrum 
also becomes soft at the SUB2 critical point A^ (see below). 

We may also solve the SUB2-m equations directly using computational techniques. In- 
deed, we study the limiting points of these approximations by using solution-tracking soft- 
ware (PITCON), which allows one to solve the SUB2-m coupled non- linear equations. We 
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FIG. 1: CCM SUB2-m and SUB2 results using the z-aligned Neel model state for the ground-state 
energy of the spin-half square-lattice XXZ model. 

again track our solution from the limit A — oo down to and beyond the limiting point 
and Fig. ^ shows our results. In particular, we note that we have two distinct branches, 
although only the upper branch is a "physical" solution. We note that the CCM does not 
necessarily always provide an upper bound on the ground-state energy, although this is often 
the case for the physical solution. By tracking from a point at which we are sure of, the 
solution we guarantee that our solution is valid, and this approach is also used for LSUBm 
approximations. We find that the two branches collapse onto the same line, namely, that of 
the full SUB2 solution, as we increase the level of SUB2-m approximation with respect to 
m. We plot the positions of the SUB2-m limiting points against l/m^ in Fig. |21 and we may 
see that these data points are found to be both highly linear and they tend to the critical 
value, Ac, for the full SUB2 equations in the limit m — > oo. It should be noted that the 
LSUBm and SUBm-m approximations also show similar branches (namely, one "physical" 
and one "unphysical" branch) which appear to converge as one increases the magnitude of 
the truncation index, m. This is a strong indication that our LSUBm and SUBm-m critical 
points are also reflections of phase transitions in the real system and that our extrapolated 
LSUBm and SUBm-m results should tend to the exact solution. 
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FIG. 2: COM SUB2-m results for the limiting points with m={ 10,40,50} for the spin-half square- 
lattice XXZ model using the z-aligned Neel model state. The unbroken line indicates a linear fit to 
these three SUB2-m limiting points, and this extrapolates to a value of 0.79854 as m — > oo. This 
result is in good agreement with the full SUB2 value for the critical point of 0.7985, indicated by 
the red star on the x = line. 

D. High-order CCM calculations using a computational approach 

We now consider the localised LSUBm and SUBm-m approximation schemes evaluated 
at high orders of approximation. We note however that lower orders of approximation may 
be evaluated by hand, as demonstrated above. However, for higher orders of approxima- 
tion we quickly find that this task becomes impossible and so we must use computational 
techniques in order to determine and solve the CCM equations. There are two methods of 
doing this. Firstly, we may determine the similarity transformed versions of the individual 
spin operators at a given level of approximation computationally and then substitute them 
into the similarity transformed version of the Hamiltonian using computational algebra, for 
example. We may however have to perform further commutations of the spin operators. 
This approach has the advantage of flexibility and we may consider any Hamiltonian, in 
principle, once the basic spin operators have been defined after similarity transformation at 
a given level of approximation. 
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Another approach is to cast the CCM ket-state correlation operator into a form given by 

TV 

s = j2s,,st + Y.Sn,A< + --- (34) 

with respect to a model state in which all spins point in the downwards 2;-direction. Note 
that Si-^^...^i^ are the CCM ket-state correlation coefficients. We now define new operators 
given by 

Gk,m = X! X! ^(^ ^ l)'5fc,m,j3, ■■■ (35) 
l>l is, 

for the spin-half quantum spin systems. We require additional terms for s > 1/2. For 
example, we may perform a similarity transform of the spin operators, which for the spin- 
half system are given by 

— 

~sl = 4 + Fkst 

= - 2F^ - iF,)ht (36) 

We now substitute these expressions into the (similarity transformed) Hamiltonian and 
evaluate the commutations - but this time by hand. The Hamiltonian is then written in 
terms of these new operators of Eq.ip^. which are themselves made up purely of spin-raising 
operators. 

We see that this approach requires more direct effort in setting up the Hamiltonian in 
terms of these new operators, whereas before we have used computer algebraic techniques 
in order to take care of this aspect. However, once this is accomplished, the problem of 
finding the ket-state equations then reduces to pattern matching of our target fundamental 
configurations to those terms in the Hamiltonian. We note that this form is then perfectly 
suited to a computational implementation because no further commutations or re-ordering 
of terms in the Hamiltonian is necessary. We note that the bra-state equations may be 
directly determined once the ground-state energy and CCM ket-state equations have been 
determined 

Q, HQ- 
Results for the ground-state energy of the spin-half square-lattice XXZ model are shown 

in Fig. Eland for the spin-half Heisenberg model (A = 1) in Table lU We note that good 
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FIG. 3: CCM LSUBm results using the z-aligned and planar Neel model states for the ground- 
state energy of the spin-half square-lattice XXZ model compared to quantum Monte Carlo results 
of Ref. Results for the LSUB6 approximation using both model states end at their respective 
critical points. 

correspondence with the results of the best and most accurate of other approximate methods 
is observed. We see clearly from Fig. 121 that the results based on the 2;-aligned model give 
the best results in the region A > 1, and those based on the planar model state appear to 
work best in the region — 1 < A < 1. This is in agreement with our understanding of this 
system that states that we have Neel ordering in the xy-plane for — 1 < A < 1 and Neel 
ordering in the z-direction for A > 1. We note that the results for the two model states are 
identical at A = 1, as we expect. 



E. Excitation Spectra of the spin-half square-lattice XXZ Model for the z-aligned 
model state 

We now consider the excitation spectrum, which may be determined using the CCM for 
the SUB2 approximation in the ground-state by assuming 

X = J2a,st . (37) 

i 

We note that the substitution of the expressions in Eqs. and (|H7|) for the ground- and 
excited-state operators, respectively, leads to an expression for the excited-state correlation 
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FIG. 4: COM LSUBm results using the z-aligned and planar Neel model states for the sublattice 
magnetisation of the spin-half square-lattice XXZ model. 



coefficients, given by 



(38) 



(We have z = 4 here.) This equation may be solved by Fourier transform techniques such 
that the an expression for the excitation spectra is given by 

1 f—:r^ ^3g^ 



e{q) = -zKJl - k^-f^{q) 



We note that K = A + 2bi, A;^ = (1 + 2A6i + 2bl)/K'^, and 7(g) = 0.5(cosqx + cosqy), and 
that we solve the SUB2 ket-state equations (for a given value of the anisotropy parameter, 
A, here) in order to obtain a value for bi. This is substituted into Eq. (jH^ above and so 
we obtain values for the excitation spectra as a function of the wave vector. We note that 
the excitation spectrum becomes soft at the CCM SUB2 critical point and results for the 
spectrum are presented in Fig. We see that the CCM excitation spectrum is identical 
in shape to those results of linear spin-wave theory (SWT) with a multiplicative factor of 
1.1672. This agrees well with results of series expansions 51 1 and quantum Monte Carlo 5^ 
that also predict a curve identical to SWT with multiplicative factors 1.18 and 1.21±0.03, 
respectively. 

Finally, we again note that the excitation energy may be determined directly from Eq. 
p2|) in "real space" without recourse to Fourier transform methods, although computational 
techniques are again necessary for higher orders of approximation. For the sake of consis- 
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3 4 5 6 1 2 

k k 

X X 

FIG. 5: Excitation spectrum for the Heisenberg model determined at the critical point at Ac = 
0.7985 for the CCM results and at A = 1 for the spin-wave theory results. The diagram plotted on 
the left is for kx = ky and the diagram on the right is for ky = 0. Note that the CCM excitation 
spectrum is identical in shape to those results of SWT with a multiplicative factor of 1.1672. This 
agrees well with results of series expansions jsil and quantum Monte Carlo ^ that also predict 
a curve identical to SWT with multiplicative factors 1.18 and 1.21ib0.03, respectively. Our results 
are in thus in good agreement with the best of other methods and this is further evidence that the 
CCM SUB2 critical point is an indication of the quantum phase transition at A = 1 in the "real" 
system. 



tency, we often retain the same level of localised approximation for the ground and excited 
states. For example, high-order CCM results are presented for the XXZ model in Fig. IHl 
and for the Heisenberg model in Table HJ We see that the CCM results converge rapidly 
with LSUBm approximation level. Indeed, extrapolated results predict that the excitation 
is gap less at A = 1, as is believed to occur for the Heisenberg model from the results of 
other approximate calculations (as discussed above). 



F. s > 1/2 Results 

We note that results for the case s > 1/2 may be generated using the CCM in the same 
way as for s = 1/2 for Neel-like model states. However, in this case we relax the condition 
that = 0, which is true for the spin-half system. We are now able to treat the 

LSUBm, SUBm-m, and SUB2 approximations in directly equivalent ways as to the spin- 
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FIG. 6: COM LSUBm results using the z-aligned Neel model state for the excited-state energy of 
the spin-half square-lattice XXZ model compared to linear spin-wave theory of Ref. [3] • 

half case. We note however that new terms are generated for the high-order formahsm and 
that this must be taken into account for s > 1/2 systems. How this is achieved is explained 

n 

in Ref. |4^. Results for the ground-state properties of the spin-one antiferromagnet on the 
square-lattice and the spin-half /spin-one (anti-)ferrimagnet are presented in Tables ITTl and 
mil respectively. We see that good results are obtained in comparison to other approxi- 
mate methods. Results using the CCM, spin-wave theory, and exact diagonalisations have 
been found to be very useful in determining the phase diagrams of the spin-half /spin-one 
ferrimagnet, even in the presence of strong frustration Clearly, there are many other 

quantum magnetic systems that can be treated in a similar fashion. 



IV. CONCLUSIONS 



The underlying formalism and practical application of the CCM have been discussed 
here. A detailed exposition of the application of the CCM to the spin-half square-lattice 
XXZ model was given, including details of analytical LSUB2 and SUB2 calculations. The 
details of applying the CCM to high orders of approximation for a localised approximation 
scheme using a computational approach were explored. Our results were seen to be in 
excellent agreement with the best of other approximate methods. The application of the 
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TABLE I: CCM results |3^, |39| for the isotropic (A = 1) spin-half square-lattice Heisenberg 
antiferromagnet compared to results of other methods. The numbers of fundamental configurations 
in the ground-state and excited-state CCM wave functions for the z-aligned Neel model state are 
given by and -^j^, respectively, and the number of fundamental configurations in the ground- 
state CCM wave function for the planar Neel model state is given by N^. Results for the critical 
points of the z-aligned Neel model state are indicated by and results for the critical points of the 
planar Neel model state are indicated by A^. (Note that results for the ground-state expectation 
values for both model states are identical for the isotropic Heisenberg model at A = 1.) 



Method 




M 






^1 






A? 


LSUB2 


-0.64833 


0.841 


1.407 


1 


1 


1 






SUB2 


-0.65083 


0.827 


1.178 








0.7985 


1.204 


LSUB4 


-0.66366 


0.765 


0.852 


7 


6 


10 


0.577 


1.648 


LSUB6 


-0.66700 


0.727 


0.610 


75 


91 


131 


0.763 


1.286 


LSUB8 


-0.66817 


0.705 


0.473 


1273 


2011 


2793 


0.843 




Extrapolated CCM 


-0.6697 


0.62 


0.00 








1.03 




LSWT[9] 


-0.658 


0.606 


0.0 








1.0 




Series Expansions [llj 


-0.6693(1) 


0.614(2) 














QMC[2] 


-0.669437(5) 


0.6140(6) 















CCM using different model states and for s > 1/2 systems was also described. It was seen 
that excellent results can be obtained for such cases. We note that we may employ the 
symmetries of the lattice and Hamiltonian in order to reduce the complexity of our problem 
- for example, by reducing the number of fundamental configurations at a given level of 
approximation. 

The application of the CCM using non-Neel model states is an interesting possibility 
in the future. Such model states might be dimer-solid or plaquette-solid model states, for 
example. This would have the advantage that novel ordering of quantum spin systems could 
be considered directly via the CCM. We would expect to obtain similar accuracy as that seen 
here using our high-order computational techniques, and this might be further extended by 
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TABLE II: CCM results for the ground state of the spin-one Heisenberg antiferromagnet at A = 1 
on the square lattice using the LSUBm approximation with m = {2,4,6}. Values for the CCM 
critical points, Ac, of the anisotropic model as a function of the anisotropy. A, are also presented. 
Note that Np indicates the number of fundamental clusters at each level of approximation. 





Np 


Eg/N 


M 


Ac 


LSUB2 


2 


-2.295322 


0.909747 


0.3240 


SUB2'' 




-2.302148 


0.8961 


0.9109 


LSUB4 


30 


-2.320278 


0.869875 


0.7867 


LSUB6 


1001 


-2.325196 


0.851007 


0.8899 


LSUBoo 




-2.3292 


0.8049 


0.98 


SWT^ 




-2.3282 


0.8043 




Series Expansions'^ 




-2.3279(2) 


0.8039(4) 





''See Ref. 
''See Ref. 
'=See Ref. 



Ill 



using parallel processing. It has recently been proven |4y, |4S| that numbers of fundamental 
clusters approaching 10^ are now possible, in principle. We might also wish to use high- 
order CCM using computational techniques for non-Neel model states in order to treat the 
excitation spectra and excitation energies. This would provide a powerful viewpoint into 
such novel-ordered states in cases where other methods might typically fail, e.g., due to 
strong frustration. Lattice quantum spin systems might provide a possible useful arena 
in which to test the properties of the extended coupled cluster method (ECCM) due 
to their underlying simplicity of their formalism and yet their complexity of behaviour for 
large numbers of particles. It would also be interesting to apply the CCM at non-zero 
temperatures, and, again, spin systems are an excellent candidate for such treatment. It 
would be interesting to see if some of the insight gained into quantum phase transitions at 
zero-temperature using the CCM might also be seen in such non-zero temperature CCM 
calculations. 
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TABLE III: Results for the ground-state energy and amounts of sublattice magnetizations 
and Mb, on the spin-one and spin-half sites respectively, of the square-lattice spin-half/spin-one 
J1-J2 ferrimagnet at Ji = 1 and J2 = based on the Neel model state. Note that Np indicates 
the number of fundamental configurations at a given level of LSUBm or SUBm-m approximation. 
CCM results are compared to exact diagonalisations of finite-sized lattices. Heuristic extrapolations 
in the limit m — > 00 are performed and a rough estimate of the error (in the last significant figure 
shown) is also given. 





Nf 


Eg/N 


Ma 


Mb 


SUB2-2 


1 


-1.192582 


0.92848 


0.85695 


SUB4-4 


13 


-1.204922 


0.90781 


0.81562 


SUB6-6 


268 


-1.206271 


0.90333 


0.80667 


LSUB6 


279 


-1.206281 


0.90329 


0.80659 


Extrapolated CCM 




-1.2069(2) 


0.898(1) 


0.796(2) 
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